/*
 * Digitizer.cc
 *
 * @created Apr 11, 2013
 * @author Tim Niggemann, III Phys. Inst. A, RWTH Aachen University
 * @copyright GNU General Public License v3.0
 */

#include "digi/G4SipmDigitizer.hh"

#include <G4DigiManager.hh>
#include <G4Timer.hh>
#include <CLHEP/Units/SystemOfUnits.h>
#include <CLHEP/Random/Random.h>
#include <CLHEP/Random/Ranlux64Engine.h>
#include <CLHEP/Random/RandExponential.h>
#include <CLHEP/Random/RandPoisson.h>
#include <CLHEP/Random/RandFlat.h>
#include <CLHEP/Random/RandBinomial.h>

#include <sys/timeb.h>
#include <math.h>
#include <set>

#include "G4SipmUiMessenger.hh"
#include "digi/G4SipmCellFireController.hh"
#include "digi/G4SipmEffectiveOvervoltageCellFireController.hh"

// G4SipmDigitizer 类的构造函数，初始化 G4SipmDigitizer 对象，并将其注册到 Geant4 的数字化管理器中。
G4SipmDigitizer::G4SipmDigitizer(G4Sipm *_sipm) : G4VDigitizerModule(_sipm->getDigitizerName()), sipm(_sipm)
{
	// Add digi collection name.
	collectionName.push_back(sipm->getDigiCollectionName());
	// Register to Geant4.
	G4DigiManager::GetDMpointer()->AddNewModule(this);
}

G4SipmDigitizer::~G4SipmDigitizer()
{
	//
}
// 该函数向给定的队列中添加串扰信号。串扰信号通过二项式过程模拟，考虑四个相邻像素。
void G4SipmDigitizer::addCrosstalk(G4SipmDigi *const d, G4SipmDigiQueue *queue)
{
	G4SipmModel *model = sipm->getModel();
	// Binomial process with 4 trials (neighbors) for crosstalk.
	// 设置邻居数量，默认为4个相邻像素。如果用户设置了不同的邻居数（例如，通过宏定义），则使用该值。
	const unsigned int nTrials = G4SipmUiMessenger::getInstance()->getNoiseCrosstalkNeighbours();
	// 计算影响范围，默认为1个像素。如果邻居数大于4，则根据公式计算影响范围。
	int shells = nTrials > 4 ? (sqrt(1 + nTrials) - 1) / 2 : 1;
	// 计算串扰概率，根据恢复校正后的概率和串扰概率，计算每个试验发生串扰的概率。
	double epsilon = 1. - pow(1. - getRecoveryCorrectedProbability(d, model->getCrossTalkProbability()), 1. / nTrials);
	if (epsilon <= 1.e-9)
	{
		return;
	}
	// 生成串扰事件数，根据二项式分布生成串扰事件的数目。
	unsigned int nCt = CLHEP::RandBinomial::shoot(nTrials, epsilon);
	// 为每个串扰事件生成新的数字化对象，并将其添加到队列中。
	for (unsigned int i = 0; i < nCt; i++)
	{
		// 生成新的数字化对象，设置类型为CROSSTALK。
		G4SipmDigi *dct = new G4SipmDigi(*d);
		dct->setType(CROSSTALK);
		CLHEP::Hep2Vector pos = model->getCellPosition(d->getCellId());
		// Dice new id.
		G4SipmCellId cellId = d->getCellId();
		// 计算新串扰像素的位置，通过随机选择x和y方向的偏移量来计算新位置。如果邻居数为4，则只允许在x或y方向上移动一个单元格。
		double a = 0, b = 0;
		while (a == 0 && b == 0)
		{
			a = CLHEP::RandFlat::shootInt(-shells, shells + 1);
			b = CLHEP::RandFlat::shootInt(-shells, shells + 1);
			// Only for four nearest neighbors.
			if (nTrials == 4)
			{
				if (a != 0)
				{
					b = 0;
				}
				else if (b != 0)
				{
					a = 0;
				}
			}
		}
		double x = a * model->getCellPitch() + pos.x();
		double y = b * model->getCellPitch() + pos.y();
		cellId = model->getCellId(x, y, false);
		if (model->isValidCellId(cellId))
		{
			dct->setCellId(cellId);
			queue->push(dct);
		}
	}
}
// 向队列中添加后脉冲事件的方法，该方法根据模型提供的概率和时间常数生成相应的后脉冲事件。
void G4SipmDigitizer::addAfterpulses(G4SipmDigi *const d, G4SipmDigiQueue *queue)
{
	G4SipmModel *model = sipm->getModel();
	addAfterpulses(d, queue, model->getApProbLong(), model->getApTauLong());
	addAfterpulses(d, queue, model->getApProbShort(), model->getApTauShort());
}
// 添加长尾后脉冲和短尾后脉冲的方法，该方法根据给定的概率和时间常数生成相应的后脉冲事件。
void G4SipmDigitizer::addAfterpulses(G4SipmDigi *const d, G4SipmDigiQueue *queue, const double probability,
									 const double tau)
{
	// Poisson distributed number.
	// 计算泊松分布的平均后脉冲数mean
	const double mean = -log(1. - getRecoveryCorrectedProbability(d, probability));
	// 根据泊松分布生成实际的后脉冲数nAp
	unsigned int nAp = CLHEP::RandPoisson::shoot(mean);
	// Bernoulli process.
	// 添加到后脉冲队列中，每个后脉冲事件的时间是通过指数分布生成的。
	for (unsigned int i = 0; i < nAp; i++)
	{
		G4SipmDigi *dap = new G4SipmDigi(*d);
		dap->setType(AFTERPULSE);
		// Dice time.
		dap->setTime(d->getTime() + CLHEP::RandExponential::shoot(tau));
		queue->push(dap);
	}
}
// 获取校正后的概率的方法，该方法根据给定的概率和权重计算校正后的概率。
double G4SipmDigitizer::getRecoveryCorrectedProbability(G4SipmDigi *const d, double prob) const
{
	// The weight equals the recharge status of the cell.
	// Measurements confirm a quadratic dependence of the probabilities to the overvoltage.
	// 击中事件的权重等于像素的恢复状态。测量表明概率与过电压的平方成正比。
	// 校正后的概率 = 原始概率 * (权重^2)
	return prob * d->getWeight() * d->getWeight();
}
// 将一系列的击中事件转换为数字化事件，并将它们添加到指定的队列中，同时记录这些事件的时间范围（最小时间和最大时间）。
void G4SipmDigitizer::addHits(G4SipmDigiQueue *queue, double *tMin, double *tMax)
{
	const G4SipmHitsCollection *hits = getHitCollection();
	if (!hits || hits->GetSize() == 0)
	{
		return;
	}
	*tMin = reinterpret_cast<G4SipmHit *>(hits->GetHit(0))->getTime();
	*tMax = *tMin;
	for (size_t i = 0; i < hits->GetSize(); i++)
	{
		G4SipmHit *hit = reinterpret_cast<G4SipmHit *>(hits->GetHit(i));
		G4SipmDigi *d = new G4SipmDigi;
		d->setType(PHOTON);
		d->setSipmId(sipm->getId());
		d->setCellId(sipm->getModel()->getCellId(hit->getPosition().x(), hit->getPosition().y(), false));
		d->setTime(hit->getTime());
		d->setWeight(hit->getWeight());
		queue->push(d);
		// Remember minimum and maximum time.
		// 记录最小和最大时间，以便后续处理。
		*tMin = std::min(*tMin, hit->getTime());
		*tMax = std::max(*tMax, hit->getTime());
	}
}
// 添加热噪声的方法，该方法根据给定的时间范围和热噪声率生成热噪声事件。
void G4SipmDigitizer::addThermalNoise(G4SipmDigiQueue *queue, const double tMin, const double tMax)
{
	// 开始计时
	G4Timer timer;
	timer.Start();
	// 检查SiPM模型的热噪声率是否为0
	if (sipm->getModel()->getThermalNoiseRate() == 0)
	{
		std::cerr << "G4SipmDigitizer::addThermalNoise(): thermal noise rate zero -> skipping." << std::endl;
		return;
	}
	// Number of cells.
	// 获取SiPM模型的像素数量
	const unsigned int nCells = sipm->getModel()->getNumberOfCells();
	// Thermal noise rate on single cell basis.
	// 计算单个单元的热噪声时间常数，热噪声的时间常数基于热噪声率的倒数计算得出
	const double tau = 1. / sipm->getModel()->getThermalNoiseRate();
	if (tau > 1. * CLHEP::second)
	{
		std::cerr
			<< "G4SipmDigitizer::addThermalNoise(): thermal noise time constant large than one second -> skipping."
			<< std::endl;
		return;
	}
	// Dice ...
	// 计算热噪声事件的起始和结束时间范围，其中包括了在噪声前和噪声后的时间。
	double tFirst = -((double)G4SipmUiMessenger::getInstance()->getNoisePreThermal() * (double)nCells) / sipm->getModel()->getThermalNoiseRate() - tMin;
	double tLast = ((double)G4SipmUiMessenger::getInstance()->getNoisePostThermal() * (double)nCells) / sipm->getModel()->getThermalNoiseRate() + tMax;

	// 生成热噪声事件，每个事件的生成时间是通过指数分布生成的。
	for (double t = tFirst; t < tLast; t += CLHEP::RandExponential::shoot(tau))
	{
		G4SipmDigi *d = new G4SipmDigi();
		d->setSipmId(sipm->getId());
		d->setCellId(CLHEP::RandFlat::shootInt(nCells));
		d->setType(THERMAL);
		d->setTime(t);
		queue->push(d);
	}
	// 输出计时信息，如果详细级别大于0
	if (verboseLevel > 0)
	{
		timer.Stop();
		std::cout << "G4SipmDigitizer::addThermalNoise(): \"" << GetName() << "\": generated " << queue->size()
				  << " thermal noise digis in (" << tFirst / CLHEP::microsecond << " us, " << tLast / CLHEP::microsecond
				  << " us) (" << timer << ")." << std::endl;
	}
}
// 模拟SiPM的击中、噪声（包括热噪声、串扰和后脉冲）以及数字化过程，最终生成并存储了数字化的SiPM信号集合。
void G4SipmDigitizer::Digitize()
{
	if (!sipm)
	{
		return;
	}
	// Start a timer.
	// 开始计时
	G4Timer timer;
	timer.Start();
	// Set verbosity.
	// 设置详细级别，根据用户界面消息传递器的实例获取的详细级别来设置。
	SetVerboseLevel(G4SipmUiMessenger::getInstance()->getVerbose());
	// Check if digitization is enabled.
	// 检查是否启用了数字化，如果未启用则直接返回。
	if (!G4SipmUiMessenger::getInstance()->isDigitizeHits())
	{
		return;
	}
	// Add signal.
	// 添加信号
	G4SipmDigiQueue queue;
	double tMin = 0, tMax = 0;
	addHits(&queue, &tMin, &tMax);
	// If queue is empty...
	// 检查队列是否为空，如果为空且未启用无信号噪声则直接返回。
	if (queue.empty() && !G4SipmUiMessenger::getInstance()->isNoiseIfNoSignal())
	{
		std::cout << "G4SipmDigitizer::Digitize(): \"" << GetName() << "\": hit collection empty." << std::endl;
		return;
	}
	// 打印击中事件的数量信息
	std::cout << "G4SipmDigitizer::Digitize(): \"" << GetName() << "\": created " << queue.size() << " photon hits."
			  << std::endl;
	// 添加热噪声，如果启用了热噪声则调用addThermalNoise方法添加热噪声事件。
	if (G4SipmUiMessenger::getInstance()->isNoiseThermal())
	{
		// Add thermal noise.
		addThermalNoise(&queue, tMin, tMax);
	}
	// Iterate queue.
	// 迭代队列，使用控制器对象来处理击中事件。根据是否模拟分流电阻来决定使用哪种控制器对象。
	double t0 = queue.empty() ? 0 : queue.top()->getTime();
	G4SipmCellFireController *controller =
		G4SipmUiMessenger::getInstance()->isSimulateShuntresistor() ? new G4SipmEffectiveOvervoltageCellFireController(sipm->getModel(), t0) : new G4SipmCellFireController(sipm->getModel(), t0);
	G4SipmDigiCollection *fired = new G4SipmDigiCollection(GetName(), GetCollectionName(0));
	// 处理队列中的每个击中事件，并根据其类型进行处理。
	while (queue.hasNext())
	{
		G4SipmDigi *d = queue.next();
		if (controller->fire(d))
		{
			fired->insert(d);
			// Handle crosstalk.
			// 处理串扰，如果启用了串扰噪声则调用addCrosstalk方法添加串扰事件。
			if (G4SipmUiMessenger::getInstance()->isNoiseCrosstalk())
			{
				addCrosstalk(d, &queue);
			}
			// Handle afterpulsing.
			// 处理后脉冲，如果启用了后脉冲噪声则调用addAfterpulses方法添加后脉冲事件。
			if (G4SipmUiMessenger::getInstance()->isNoiseAfterpulse())
			{
				addAfterpulses(d, &queue);
			}
		}
	}
	// Store digis.
	// 存储数字化事件
	StoreDigiCollection(fired);
	// Free space.
	// 释放空间
	delete controller;
	// Stop timer and print information.
	// 停止计时器并打印信息。
	timer.Stop();
	std::cout << "G4SipmDigitizer::Digitize(): added digi collection \"" << fired->GetName() << "\" for digitizer \""
			  << GetName() << "\" with " << fired->GetSize() << " entries (" << timer << ")." << std::endl;
}
// 从数字化管理器中获取一个与特定SiPM相关的命中集合。
const G4SipmHitsCollection *G4SipmDigitizer::getHitCollection()
{
	// We search and retrieve the hits collection.
	// 查找并检索击中集合。
	G4DigiManager *digiManager = G4DigiManager::GetDMpointer();
	G4int collId = digiManager->GetHitsCollectionID(sipm->getHitCollectionName());
	if (collId >= 0)
	{
		const G4SipmHitsCollection *hitsCollection =
			static_cast<const G4SipmHitsCollection *>(digiManager->GetHitsCollection(collId));
		return hitsCollection;
	}
	return NULL;
}
// 获取SIPM对象的方法，该方法返回与数字转换器关联的SIPM对象的指针。
G4Sipm *G4SipmDigitizer::getSipm() const
{
	return sipm;
}
// 设置SIPM对象的方法，该方法将给定的SIPM对象的指针赋值给sipm成员变量。
void G4SipmDigitizer::setSipm(G4Sipm *_sipm)
{
	sipm = _sipm;
}
